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Motivated by the current interest in the quantum dimer model on the triangular lattice, we 
investigate the phase diagram of the closely related fully-frustrated transverse field Ising model on 
the honeycomb lattice using classical and semi-classical approximations. We show that, in addition 
to the fully polarized phase at large field, the classical model possesses a multitude of phases that 
break the translational symmetry which, in the dimer language, correspond to a plaquette phase 
and a columnar phase separated by an infinite cascade of mixed phases. The modification of the 
phase diagram by quantum fluctuations has been investigated in the context of linear spin-wave 
theory. The extrapolation of the semiclassical energies suggests that the plaquette phase extends 
down to zero field for spin 1/2, in agreement with the \/T2 x \/l2 phase of the quantum dimer model 
on the triangular lattice with only kinetic energy. 
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I. INTRODUCTION 

Quantum dimer models have emerged as one of the 
main paradigms in the investigation of quantum spin liq- 
uids. The Rokhsar-Kivelson (RK) quantum dimer model 
(QDM), which includes a potential interaction of ampli- 
tude V between dimers facing each other and a kinetic 
term of amplitude t flipping them around rhombic pla- 
quettes, has recently attracted special attention. The 
main reason comes from the presence on the triangular 
lattice of a resonating valence bond (RVB) phase first dis- 
covered by Moessner and Sondhi^ and extensively studied 
since then using zero temperature Green's function quan- 
tum Monte Carlo (GFQMC).^"'* Exact results have been 
obtained at the RK point {V/t — 1), where the sum of all 
configurations can be proven to be a ground state, ^ and 
at V > t, where the non-flippable configurations are the 
ground states. Analytical results have also been obtained 
in the limit V/t — > — oo, where columnar states have been 
shown to be selected. However, in the intermediate range 
below the RK point, most of what is known about the 
model is based on numerical simulations. 

A closely related model for which a number of ana- 
lytical results have already been obtained is the fully 
frustrated transverse field Ising model (FFTFIM) on the 
honeycomb lattice defined by the Hamiltonian: 

where F > is the transverse magnetic field, J > 
is the coupling constant of the Ising interaction term, 
denotes pairs of nearest neighbors on the honey- 
comb lattice, and Mij = ±1 is such that for each hexagon 
of the lattice the number of antiferromagnetic bonds 
{Mij — —1) is odd, different choices of My correspond- 
ing to the same model up to the rotation of some spins 
by TT around the x axis^. Transverse field Ising models 



have been the subject of intense investigations over the 
years. ^ The relationship between the FFTFIM on a reg- 
ular lattice and the QDM on the dual lattice has been 
first emphasized by Moessner, Sondhi and Chandra^ who 
showed (see also Rcf. 1) that, in the limit F/J — >■ 0, the 
FFTFIM on the honeycomb lattice maps onto the QDM 
on the triangular lattice with t = F^/J and V = Q. For 
the FFTFIM on the honeycomb lattice, they also car- 
ried out a Landau-Ginzburg analysis and identified four 
soft modes which, upon lowering F/J, simultaneously be- 
come gapless, leading to a surprisingly large unit cell of 48 
sites. Details of this calculation have been reported later 
by Moessner and Sondhi in Ref. 8. These authors fur- 
ther conjectured that the translational symmetry break- 
ing transition out of the paramagnetic phase coming from 
large F/J provides a reasonable description of the transi- 
tion between the RVB phase and the intermediate phase 
of the QDM on the triangular lattice.^ 

Building on this conjecture, Misguich and one of the 
present authors have carried out a semiclassical investi- 
gation of the paramagnetic phase of the FFTFIM on the 
honeycomb lattice*^ and have shown that the dispersion of 
the spin waves and their softening at the transition are 
in remarkable agreement with the dispersion of visons 
in the QDM on the triangular lattice and their crystal- 
lization transition as revealed by quantum Monte Carlo 
(QMC) simulations."* However, the analysis of Ref. 9 has 
not covered the small F/J parameter range. 

In the present paper, we perform a systematic investi- 
gation of the FFTFIM on the honeycomb lattice in the 
complete parameter range < F/J < -|-oo with classical 
and semi-classical approximations. As we shall see, the 
classical phase diagram is much richer than expected, 
with an infinite number of different crystalline phases 
below the paramagnetic phase: a plaquette phase, a cas- 
cade of mixed phases, and a highly degenerate columnar 
phase. Quantum fluctuations have been treated within 
linear spin-wave theory, leading to a partial lifting of the 
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FIG. 1: Sketch of the gauge used in most of the paper. Here 
and below antiferromagnetic bonds (with Mij — —1) are 
shown by zigzags, all other bonds being ferromagnetic (with 
Mij = +1). 

degeneracy of the columnar phase, and to an increase of 
the size of the region occupied by the plaquette phase. 

To make contact between the physics of the FFTFIM 
and of the QDM, it is useful to introduce a gauge theory 
defined on the triangular lattice by the Hamiltonian: 

H = -JT.-t-r Y^Ii-m^ (2) 

( i Hi) 

where i runs over the sites of the dual honeycomb lattice, 
and are the three bonds forming the triangular pla- 
quette around site i. As shown by Moessner, Sondhi and 
Fradkin,^" the FFTFIM is equivalent, up to a twofold de- 
generacy, to the odd sector of this gauge theory defined 

by 

/[a] 

for all a, where a is a site of the triangular lattice, and 
the product over l[a] runs over the six links emanating 
from a. For a compact discussion of the correspondence 
between the three models, see e.g. the introduction of 
Ref. 9. 

The discussion of the ordered phases is simpler in the 
context of the gauge theory. Indeed, in the FFTFIM lan- 
guage, the actual orientation of the spins in a given state 
depends on the choice of the matrix My . By contrast, 
the dimcr operator of the gauge theory defined by 

di^^il-rn, (4) 

translates into 

in the Ising language, and its expectation value does not 
depend on the choice of Mij. Another advantage of the 
gauge-invariant language is that it allows to make a di- 
rect comparison with the numerical results obtained on 
the QDM since they live on the same lattice and are 
defined in terms of the same link operators. So, while 
all reasonings and calculations will be performed in the 
context of the FFTFIM, the only formulation adapted 
to the semiclassical approach, the structures of different 



ordered phases will be also described in gauge-invariant 
terms. Throughout the paper, we use only gauges in 
which each hexagon of the lattice contains exactly one 
antiferromagnetic bond (with Mij = — 1) and five ferro- 
magnetic bonds (with Alij — +1). Most results will be 
presented for the simplest periodic arrangement of the 
antiferromagnetic bonds shown in Fig.l. However, this 
choice of gauge does not always lead to the smallest possi- 
ble unit cell in terms of the spin representation. Thus, we 
will also introduce other gauges whenever this is helpful. 

The paper is organized as follows. In section II, we 
concentrate on the limit T/J <C 1, which has not been 
considered in Refs. 1,7-9, and we show that columnar 
phases reminiscent of the V —?' —oo limit of the QDM 
are stabilized. In section III, we revisit the vicinity of the 
RVB phase. We recover the symmetry predicted by the 
Landau-Ginzburg approach of Refs. 7,8 and by the spin- 
wave analysis of Ref. 9, but we find that the bonds with 
the largest dimer density form separate 4-site rhombic 
plaquettes instead of having a uniform distribution inside 
a 12-site unit cell as reported in Ref. 9. The reasons for 
this discrepancy are explained in subsection IIIC. In sec- 
tion IV, we discuss the transition between the plaquette 
phase and the columnar phase and show that they are 
separated by a region of intermediate phases of mixed 
character. The stability of these phases with respect to 
quantum fluctuations and the semi-classical phase dia- 
gram are discussed in section V. The paper ends with a 
short conclusion in section VI. 



II. COLUMNAR PHASE 

In this section we discuss the properties of the model 
when T/J is small. The argument proceeds in three 
steps. First we determine the ground state manifold of 
the Heisenberg model with purely Ising-like interactions 
in the absence of magnetic field (F = 0). Then we in- 
vestigate how the extensive degeneracy of these ground 
states is lifted by a small transverse field. Finally, we 
discuss the effect of quantum fluctuations in the context 
of linear spin- wave theory. 

A. Zero transverse field 

In the absence of a transverse magnetic field (F = 0) , 
we are left with a model without quantum fluctuations 
in which the interaction term couples only the z com- 
ponents of neighboring spins on the honeycomb lattice. 
With our choice of gauge, one bond on each hexagon is 
antiferromagnetic {Mij = — 1) and the others are ferro- 
magnetic {Mij = 1). Frustration is present since it is 
clearly impossible to minimize the energy of all bonds of 
a given hexagon. 

For Ising spins, i.e. spins which can only point up or 
down along the z direction, the best one can do is to 
satisfy five bonds leaving one unsatisfied. This can be 



3 



done in six different ways according to which bond is 
not satisfied ("frustrated") and the resulting energy is 
—4 J. Up to a global reversal of the spins, a ground state 
is characterized by the distribution of frustrated bonds 
such that there is exactly one of them per hexagon. 

For three-dimensional vectors of norm S, the situation 
is slightly more subtle because the twelve Ising configu- 
rations with all spins parallel or antiparallel to z axis are 
not the only ground states of a single hexagon. To see 
this, let us consider a single hexagon and investigate the 
possibility of a given spin i not to be directed along z. 
The variation of the energy of the hexagon 



B. Classical ground states in small transverse field 

Let us now switch on a small transverse field and study 
how the local minima of the classical Hamiltonian evolve 
upon increasing the field. Since the field is along x, 
the spins are expected to acquire a small x component, 
and to describe the spin configuration evolving from a 
given ground state of the pure Ising case, we use the 
parametrization: 



Sf = S'sin6',, 
Sf = (7,5 cos ( 



(9) 



hex 



i=i 



1 cos 9j cos 9 



(6) 



(where the angle 9j parameterizes the deviation of spin j 
from the z axis) with respect to 9i leads to the condition 



Mi_i,i cos 6ij_i + Ali^i+i cos 6ij+i = . 



(7) 



If this condition is satisfied, the terms in Eq. (6) which 
depend on 9i drop out, so that one is left with the energy 
of an open chain of five spins. In an open chain one can 
trivially minimize the energy of each bond by choosing 
cos0j_i_i = Afj.j+i cos 6'j = ±1 which leads to = — 4J 
and to the automatic fulfillment of condition (7), leaving 
9i arbitrary. Note that this argument excludes a devi- 
ation from the z axis of more than one spin, since the 
energy of a five-spin open chain cannot be as low as —4 J 
if not all five spins are along z. So, for three-dimensional 
spins, the energy of a single hexagon is minimal as soon 
as it is minimal for four consecutive bonds, and the spin 
at the remaining site can have any direction. 

It is natural to ask whether this additional freedom 
increases the degeneracy of the ground state manifold 
of the continuous model in comparison with the case of 
Ising spins. To demonstrate that this is not the case, let 
us assume that at site i the spin is not along z. To mini- 
mize simultaneously the energy of the three hexagons to 
which it belongs, three conditions of the form (7) must 
be fulfilled: 



=0, 

1^2 + ^3 = 0, 

Y3 + Yi=0, 



(8) 



where Yq = M^ ^^ cos9i^ = ±1 (with a = 1,2,3) and ia 
are the three nearest neighbors of site i. It is evident 
that the restriction Ya — ±1 does not allow all three 
equations (8) to be satisfied simultaneously. Therefore, 
it is impossible for any spin not to point along z, and the 
ground state manifold coincides with that of the frus- 
trated Ising model with the same lattice, i.e. it consists 
of all Ising configurations with one frustrated bond per 
hexagon. Each of these states is a local minimum of the 
Hamiltonian. 



where ct^ = ±1 is the sign of S'f and is determined by 
the ground state of the pure Ising case around which we 
expand. In terms of the gauge-invariant bond variable 
Tij = MijaiUj, which is equal to -1(+1) if the bond 
is frustrated (not frustrated) , the classical energy can be 
rewritten as 



E ■ 



Ti, COS Ui COS ( 



r^sin( 



(10) 



In the limit F ^ J the deviations from the z di- 
rection are small, and the classical energy can be ex- 
panded in the variables 9i around Oi = 0. To sec- 
ond order, the interaction term in Eq. (10) decouples: 
cos 01 cos 0j K. Tij{\ — 9l/2 — 6'|/2). Now, for any 
ground state of the pure Ising case, the set {Tij} is such 
that only one bond in each hexagon is frustrated. There- 
fore each site belongs at most to one frustrated bond. If 
we denote by F (resp. NF) the set of what we call below 
frustrated (resp. non frustrated) sites, namely, the sites 
belonging to one frustrated bond (resp. no frustrated 
bond), the energy up to second order can be rewritten: 



= i?r=o + E Ue^ - T9^ + ^ ( - ve, 



ieNF 



Minimizing E^"^^ with respect to {9i\ leads to 



r F/j 

\F/3J 



for « e F , 
for i e NF . 



(11) 



(12) 



Since the number of frustrated and non frustrated sites 
is the same for all ground states, the energy up to second 
order in 9i is the same in all ground states. So, second 
order corrections do not lift the degeneracy. They only 
induce a difference in orientation between the spins which 
belong to a frustrated bond and those which do not. 

So to lift the degeneracy we have to push the expansion 
in 9i to higher orders. To fourth order, it reads: 



£;(4) 



i6NF L 



3J 



J 



. r\2n2 



F e. 



4! 



4! 



-V[9,- 



3! 



3! 



(13) 
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From the previous discussion, we know that the values 
of di minimizing the energy to order 0(9^) are given by 
Eq. (12). Injecting these solutions into the fourth order 
expansion of the energy, we notice that the terms 9f and 
9f only contribute in two different ways depending on the 
type of site (frustrated or non frustrated) . They will thus 
not lift the degeneracy. By contrast, the crossed terms 
Tij9f9'j contribute in four different ways depending on 
the environment of the sites i and j. The four cases are 
illustrated in Fig. 2. 

The contributions of the fourth order crossed terms 
to the energy for the different configurations in units of 
rV4J3 are +1 for 2(a), -i for 2(c), for 2(b) and 
— 1 for 2(d). Since these energies are not equal, these 
crossed terms are expected to lift the degeneracy, at least 
partially. 



(a) 



(b) 



1 J \ 



(c) 



(d) 



FIG. 2: (Color online) Local configurations of frustrated 
bonds leading to different contributions of tfie fourth order 
crossed term 
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For a lattice of A'hex hexagons the total number of 
bonds is SiVhex- The constraint that each hexagon has 
one frustrated bond implies that the number of frustrated 
bonds is equal to A'hcx/2. This fixes the number Na of 
configurations 2(a) to be equal to A'^hGx/2. By contrast, 
the numbers of configurations of type 2(b), 2(c) and 2(d) 
(respectively Nh,Nc,Nd) depend on the way the frus- 
trated bonds are arranged on the lattice. However Nb, Nc 
and are not independent but have to satisfy the fol- 
lowing relations: 



Nb + N, + Nd 



hex : 



(14) 



(15) 



Eq. (14) comes from the conservation of the total number 
of bonds Na + Nf, + N^ + N^ = SiVhox, whereas the right- 
hand side of Eq. (15) comes from counting all frustrated 
bonds by looking at how many of them are adjacent to 
each of the non frustrated bonds. The result of this cal- 
culation has to be divided by four, because in such a 
procedure each frustrated bond is counted four times. 



The total contribution of the fourth order crossed 
terms of the energy can then be written as: 



ii'j)'{Na-N,- 



J (r\^ C^N — 64 A7 
Vsi 81 



(16) 

This contribution is a decreasing function of Nd, so the 
lowest energy will be reached for the largest possible 
value of Nd- Now, since there is only one frustrated 
bond per hexagon, Nd cannot exceed the number of frus- 
trated bonds, Na — iVhox/2. This upper limit is reached 
for configurations in which all the frustrated bonds are 
organized into chains of alternating frustrated and non 
frustrated bonds (see examples in Fig. 3). In what fol- 
lows we refer to this family of states as columnar states 
(see Fig. 3). In columnar states Eqs. (14) and (15) fix 
both Nb and Nc to be equal to A'^hcx- 

So, the fourth order contribution to the energy par- 
tially lifts the degeneracy and selects the family of colum- 
nar states. A priori, higher orders might further lift the 
degeneracy. That this is not the case is best seen by 
constructing the exact local minima that correspond to 
columnar states. We start by rewriting the energy: 



E 



J 



J 



cos 9i (cos 9i^ + cos 9i^ + cos i 



r sin 9j 



cos 9j (— cos 9j^ + cos 9j^ + cos 9j^ ) + F sin 9j 

(17) 

where ii,i2, is (resp. Ji, ^2, js) are the three neighbors of 
site i (resp. j), and the frustrated bond is taken to be 
between sites j and ji. To minimize the energy, the set 
of angles {6'i, 9j} must be a solution of the equations: 

^ — J sin 9i (cos -I- cos 9i^ + cos ) — F cos 6*^ = , 

^ = J sin 9j {— cos 9j^ + cos 9j^ + cos 9j^ ) — F cos 9j = . 

(18) 

Now, in columnar structures, all frustrated sites have 
identical environments (with exactly two frustrated 
neighbors) and all unfrustrated sites also have identical 
environments (with exactly one frustrated neighbor). So, 
if the angles 9i and 6*2 satisfy the equations: 



J sin 6*1 (2 cos 6*1 
J sin 02 cos 6'i — F cos 02=0 
then the set of angles 



COS6I2) - F cos 6*1 = 0, 



91 for i e NF 

92 for i e F 



(19) 



(20) 



is a solution of Eqs. (18). The non trivial solutions of 
Eqs. (19) describing the evolution of columnar states with 
the change of F/J are given by: 



smfi = 



_ sin(/?/3) 



cos(;3) 



sm 92 



_ sin(/3) 
~ cos(;3/3) 



(21) 
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where tan j3 = T / J . 

The substitution of Eq. (20) into Eq. (10) shows that 
the classical energy of a columnar state is given by 

N 

Ecoi = — — [J cos 01 (cos 6*1 + COS02) + r(sin^?i + sin 6*2)] , 

(22) 

where N is the total number of sites. Naturally, the vari- 
ation of Ecoi with respect to 9i and 62 reproduces Eqs. 
(19) which we used to find the values of 9i and ^2- In or- 
der to verify that it never becomes more advantageous to 
minimize Nd rather then to maximize it, we have studied 
also the solutions with Nj, — and checked that for any 
relation between T and J they have higher energy than 
the columnar states (see Appendix A). 

A convenient classification of columnar states can be 
introduced by describing them in terms of zero-energy 
domain walls formed on the background of the sim- 
plest columnar state, an example of which is shown in 
Fig. 3(a). Below we call it the I''* columnar state. In this 
state all frustrated bonds have the same orientation and 
form straight columns shown in the figure by the shading. 
In terms of Fig. 3 the walls of the first type are horizontal 
and take place whenever the orientation of the frustrated 
bonds changes from left to right. The 2""^ columnar state 
(Fig. 3(b)) corresponds to the configuration having the 
highest possible density of such domain walls. 

The domain walls of the second type are perpendicu- 
lar to the frustrated bonds, and correspond to changing 
the orientation not of frustrated bonds but of columns. 
The 3'''^ columnar state (Fig. 3(c)) is the configuration 
having the highest possible density of walls of the second 
type as the orientation of the columns changes at every 
frustrated bond. Other columnar states having the same 
classical energy can be obtained by introducing arbitrary 
sequences of parallel domain walls either of the first or 
of the second type separating domains of the 1*^' colum- 
nar state. An analogous classification of columnar states 
had earlier been introduced by Moessner and Sondhi^ in 
terms of the QDM. 

Fig. 4(a) presents a plot of the dimer density for the 
1**' columnar state at F/J = 1.5. The bonds of the dual 
triangular lattice having the highest dimer densities are 
organized into a columnar pattern. Fig. 4(b) is a plot of 
the I*'' columnar state in the classical spin model. The 
component of the spin on frustrated sites (green arrows 
in Fig. 4(b)) is smaller than that on non frustrated sites. 

C. Quantum fluctuations 

The effect of quantum fluctuations on the columnar 
states, in particular their local stability and their de- 
generacy, has been investigated in the context of linear 
spin- wave theory (LSWT). It is impossible to perform a 
LSWT calculation for all columnar states since the fam- 
ily is infinite and contains many members which are not 
periodic. The logic we have followed is based on the ex- 
pectation that the difference in energy between each pair 




(a) l'^' columnar state (b) 2"'' columnar state 




(c) 3'^'^ columnar state (d) 4"^ columnar state 



FIG. 3: (Color online) Examples of columnar states. The 
frustrated bonds are represented as dashed red lines. In the 
dimer representation, the bonds of the dual triangular lattice 
which intersect the frustrated bonds of the honeycomb lattice 
have the highest dimer density. The 4*'^ columnar state dilTers 
from the S"^"* one by having exactly half the number of domain 
walls of the second type (see main text). 



of states is determined primarily by the difference in the 
number of domain walls they contain. 

In Sec. II B we have established that the structure 
of columnar solutions is described by Eqs. (9), where 
<7i—±.l is determined by the ground state of the pure 
Ising case and the values of the variables Oi are given by 
Eqs. (20) and (21). It is convenient to start the construc- 
tion of the Hamiltonian describing the harmonic fluctua- 
tions around these states by performing a rotation of the 
spins on each site, 

Sf = a,; cos 6l,5f' + sin 61, S'f, 

SI = Sf, (23) 

Sf = -sine.Sf + <7iCose^SI' . 

in such a way that the Hamiltonian expressed in terms 
of the variables S^' and S'^' has a ferromagnetic ground 
state. 

Mapping the new spin operators to Holstein-Primakoff 
bosons in the harmonic limit 

St ^S- , sf « (a. + 4) , (24) 
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\/\/\/\/\/\/ 

\/\/\/\/\/\/ 

(a) l'^* Columnar state: dimer 



representation 



a 



(b) l'^' Columnar state: spin representation 

FIG. 4: (Color online) (a) Plot of the dimer density dij at 
r/J = 1.5 for the columnar state. The thickness of the 
bonds is proportional to dij . The dark blue bonds correspond- 
ing to the highest dimer density are organized into columns, 
(b) Spin configuration in the 1*^* columnar state in the gauge of 
Fig. 1 (with the same notation for antiferromagnetic bonds). 
The two types of arrows correspond to the two spin orienta- 
tions realized in that state. The unit cell is defined by the 
two vectors a and b with 161 — \'^\a\. 



then yields the quadratic Hamiltonian: 

H = -Ecoi + 7i X! "^"^ + X! "-l^-i 
J 



(25) 
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25* 



Mij sin di sin 0j ^aiUj + a|aj + h.c. 



' cos 6*2) +r sin 6*1] , (26) 
(27) 



where 

71 = (l/S')[Jcos6'i(2cos6li 

72 = (l/S')(Jcos6'iCOs6l2 + rsin6'2) , 

and E'j.qj is the classical energy of a columnar state. 
Eq. (25) can be reduced to a gauge-invariant form (with 
Mij replaced by Tij) by replacing by aiUi and a| by 
CTja^ in Eqs. (24). However, we use Eq. (25) in the follow- 
ing because it allows an easy proof that domain walls of 
the first type do not change the energy of the harmonic 
fluctuations. 

It is evident that for 6i given by Eq. (20), the expres- 
sion in the right-hand side of Eq. (25) is exactly the same 
for all columnar states having the same sets of frustrated 
and non frustrated sites. Since the introduction of do- 
main walls of the first type interchanges only the po- 
sitions of frustrated and non frustrated bonds forming 



straight columns, but does not change the positions of 
frustrated sites [see Fig. 3(a) and Fig. 3(b)], the expres- 
sion in the right-hand side of Eq. (25) will be exactly the 
same for all columnar states which can be transformed 
one into another by the introduction of some number 
of domain walls of the first type. This proves that the 
contribution of the harmonic fiuctuations to the energy 
is the same for all members of the family of columnar 
states having only domain walls of the first type. 

After partitioning the honeycomb lattice into four sub- 
lattices in accordance with the structure of the unit cell 
shown in Fig. 4(b) and performing on each sublattice the 
Fourier transformation with wavevector q, the quadratic 
bosonic Hamiltonian of the T** columnar state is reduced 
to the form 



H = E, 



col 



7V^[4F(g)a,--(7i+72)] 



(28) 



In this expression, a- is an eight-component vector 
{a-if,i,a_if^2,a-q,3,a-ff^4,ali,at.^,atj,at j, where a^,„ 
are the bosonic operators with wavevector q acting on 
the n^^ sublattice, and H{q) is an 8 x 8 hermitian matrix 
given by: 
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(29) 



where 



A*EE/i(g)--^^(-l + e^«"^), 
^^rj(q) = -l^{l + e^?S), 



(30) 



2S 



The vectors a and b are shown in Fig. 4(b). 

As discussed above, the harmonic Hamiltonian is the 
same for the whole family of columnar states constructed 
by introducing an arbitrary number of domain walls of 
the T'* type. This family includes for instance the 2"^^ 
columnar state. In the harmonic approximation, all these 
states have the same quantum corrections to the energy, 
therefore to order 1/S the degeneracy is not lifted. Note 
however that the absence of degeneracy lifting for this 
family of states at the harmonic level is not related to 
a symmetry of the original Hamiltonian. So we expect 
this degeneracy to be removed if one goes beyond the 
harmonic approximation, and higher order terms are ex- 
pected to select either the I''* or the 2"*^ columnar state 
depending on whether the energy of a domain wall of 
the first type is positive or negative. However the ef- 
fect of anharmonicities has not been investigated in this 
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work. Note that a similar effect, namely the incapac- 
ity of harmonic fluctuations to fully lift a well-developed 
accidental degeneracy of the ground states, has already 
been reported for various other models (in particular, 
with kagome,^^^^'^ honeycomb, dice^^ and pyrochlore^^ 
lattices) . 

By contrast, the 3'^'^ columnar state is described by 
a different harmonic Hamiltonian which is not written 
down here explicitly because the number of sites per unit 
cell, hence the linear dimension of the matrix H(q), is 
twice as large, so that the matrix H{q) is 16 x 16. The 
energy of zero point fluctuations in this state turns out 
to be higher than in the I''' columnar state (see Fig. 5). 
This suggests that domain walls of the second type have 
a positive energy. 

To support this statement, we have applied the same 
reasoning as used in Ref. 18 for the investigation of the 
frustrated XY model on a triangular lattice and have 
considered the 4'^ columnar state (Fig. 3(d)) which dif- 
fers from the 3'''^ one in that the density of domain walls 
of the second type is exactly half as large. Fig. 5 com- 
pares the numerically calculated differences between the 
value of the quantum corrections to the energies of the 
2nd^ 3rd ^j-^j ^th columnar states and its value for the 1*** 
columnar state. In particular, the inset in Fig. 5 presents 
the ratio of these quantities for the 3'''^ and 4**^ states. 
This ratio is very close to two, supporting the suggestion 
that the fluctuation induced corrections to the energy are 
essentially proportional to the density of domain walls of 
the second type. 



^/site 
1.5-10-^ 



1-10-3 
5-10-^ 




.^2 





« C2 —C\ , 
'C:,-Ci : 
* C4 ~ C\ • 



0.5 1 1.5 2 

r/j 



0.5 1 1.5 

r/j 



FIG. 5: (Color online) Energies (per site) of the 2""^ (red 
crosses), 3'^'^ (green circles) and 4"' (orange diamonds) colum- 
nar states calculated in the harmonic approximation, counted 
with respect to the energy of the 1*^' columnar state and ex- 
pressed in units of J. The inset is a plot of the ratio of the 
energy of the 3'^'^ columnar state over that of the 4**^ columnar 
state. 



Upon increasing F/J, the classical states remain lo- 
cally stable until soft-modes appear in the spin- wave dis- 
persion. For all columnar states without domain walls of 
the second type this takes place at T/jKi 2.004, and for 
the 3"'^ columnar state at F/J « 2.373. To summarize. 



harmonic fluctuations partially lift the degeneracy of the 
classical ground state manifold in favor of the columnar 
states having only domain walls of the F** type. 



III. PLAQUETTE PHASE 

A. Soft modes and the ground state periodicity 

In the limit J = the Hamiltonian consists simply of a 
coupling to the transverse magnetic field F, and the clas- 
sical ground state is completely polarized with all spins 
aligned along the magnetic field in the x direction. The 
same state minimizes the classical energy for sufhcicntly 
large ratio F/J. With the choice of gauge of Fig. 1, the 
unit cell of this state contains 4 sites (see Fig. 6). 




FIG. 6: (Color online) In the polarized state all spins are 
aligned along the magnetic field. The unit cell of this state 
is the same as that of the 1"* columnar state: it is defined by 
the vectors a and b. 

The analysis of Refs. 7-9 indicates that the polar- 
ized phase becomes unstable at F = Fc = V6J. At 
this value of the field, soft modes appear in the disper- 
sion relation at momenta {qx,(lz) — ^ ( gjoi ' 2jl)| ) ^'^^ 
{qx,qz) = ±(^: 2|E|)' triggering a second order transi- 
tion to a new phase whose periodicity can be determined 
from the q points corresponding to the soft modes. 

Since any linear combination of these four modes is 
invariant under translations by vectors (3a — b) and 46, 
a state associated with them should have the periodic- 
ity in real space imposed by these two vectors that de- 
fine a unit cell containing 48 sites of the honeycomb lat- 
tice (Fig. 7(b)). Moreover, since any linear combination 
of the four soft modes under the translation by 26 just 
changes sign, this cell should allow a division into two 
halves which in the spin representation differ from each 
other only by the reflection of all spins about the x axis 
but in terms of gauge-invariant variables are identical. 

There exists a possibility to make these two halves re- 
ally identical in terms of spin representation as well just 
by choosing a different gauge shown in Fig. 7(c). In this 
gauge a state related to the soft modes listed above is 
periodic with a 24-site unit cell defined, for example, by 
vectors 3a — 6 and 26. However, if one uses the simplest 
gauge of Fig. 1 and imposes periodic boundary conditions 
along the x and z directions, the periodicity dictated by 
the wave vectors of the soft modes requires to use a cell 
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of size 12a x 4& that contains 192 sites of the honeycomb 
lattice.* 



B. Numerical minimization of energy 

The minimization of the classical energy using Math- 
ematica minimization routines for the 192-site system 
with periodic boundary conditions have confirmed that 
the real periodicity of the classical ground state in the 
gauge of Fig. 1 is determined by a 48-site unit cell which 
can be divided into two halves in such a way that the 
second half differs from the first one by the reflection of 
all spins about the x axis. Inside the cell one finds a 
pattern of six different orientations of the spins as well 
as their reflections about the direction of the field. 

The structure of the state minimizing the classical en- 
ergy is shown in figure 7(b). The radii of the circles are 
proportional to the absolute value of the z component of 
the spins 15^1 and the different signs of S*^ are kept track 
of by plotting full and empty circles. is not plotted 
but is always positive since the spins tend to align with 
the magnetic field. The size of the elementary cell can 
be reduced to 24 sites by choosing the gauge depicted in 
Fig. 7(c) by zigzagged bonds. In this gauge the sign of 
is the same for all spins and the spin pattern is centered 
on one of the sites of the honeycomb lattice. 

Naturally, it is even more convenient to discuss the 
structure of an ordered state in terms of gauge-invariant 
dimer density dij defined by Eq. (5). In the polarized 
phase (at T/J> -\/6), Sf = for all sites i, so that the 
dimer density is uniform and equal to ^ on all bonds. 
Below the critical magnetic field, = VQJ, the dimer 
density on many bonds becomes smaller than ^ . For the 
pattern of dij the two halves of the 48-site elementary cell 
are identical because the dimer pattern is conserved when 
reversing the sign of for all spins. Accordingly, the 
elementary cell corresponds to 24 sites of the honeycomb 
lattice or to 12 sites of the triangular lattice dual to it. In 
other terms, the periodicity of the dimer density pattern 
is the same as in the a/T2 x \/12 phase found around 
V/t = in the QDM on the triangular lattice. 

In Fig. 7(a) the elementary cells are represented by the 
large hexagons. Since inside an elementary cell the dimer 
density plot displays a pattern of four-site plaquettes hav- 
ing the highest dimer density (see Fig. 7(a)), following 
the convention adopted in the QDM literature^^ we refer 
to this phase as the plaquette phase. This phase is the 
analog of the a/12 x -^12 phase found around V/t = in 
the QDM.i 

Note that the dimer density plot obtained below 
r/j = V6 in our calculation (Fig. 7(a)) differs signifi- 
cantly from the one presented in Ref. 9. The two plots 
have the same symmetry, P31m, but the pattern of Ref. 9 
does not reveal four-site plaquettes. In fact, the differ- 
ence can be traced back to the fact that the solution of 
Ref. 9 was obtained by a variational calculation in the 
subspace of linear combinations of the four soft modes 




(a)Plaquette phase: dimer (b)Plaquotto phase: unit 
representation cell in the spin 

representation 




(c)Plaquette phase: the smallest 
unit cell in the spin representation 



FIG. 7: (Color online) (a) The dimer density dij in the plaque- 
tte phase at r/J = 2. The thick blue bonds corresponding to 
the highest density (dij = |) are organized into 4-site rhom- 
bic plaquettes. On all other bonds the dimer density satisfies 
< dij < ^ , the thickness of the bonds being proportional to 
dij . 

(b) The spin configuration in the same state in the gauge of 
Fig. 1. The radii of the circles are proportional to \Si\, while 
positive and negative values of Si are represented as full and 
empty circles. The green dashed rectangle shows the 48-site 
unit cell (3a — b) x 4fo. It can be split into two halves which 
differ from each other by the sign of 5^. 

(c) The same spin configuration in the gauge that leads to a 
24-site unit cell (large hexagon). As before antiferromagnetic 
bonds fixing the gauge are depicted as zigzag bonds. The 
sites at which the classical spins have the same values of Sf 
are labelled with the same number. Note the existence of 6 
sites with Sf = 0. 

(which minimize the sum of the second and fourth or- 
der contributions to the classical energy), whereas the 
present solution was obtained by assuming that the soft 
modes dictate only its periodicity. The reason why the 
two solutions do not have the same asymptotic form when 
r/J tends to \/6 from below is detailed below in Sec. 
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IIIC devoted to the analytical investigation of the pla- 
quette state structure in the vicinity of the phase transi- 
tion. 

The degeneracy of the plaquette phase is equal to 48 in 
terms of the spin representation and to 24 in terms of the 
dimer representation. Each of the 24 equivalent dimer 
patterns [one of which is shown in Fig. 7(a)] corresponds 
to two spin configurations which can be transformed one 
into the other by changing the sign of S'f for all spins. 

The local stability of the plaquette phase with respect 
to quantum fluctuations has been investigated within the 
gauge of Fig. 7(c) to reduce the hermitian matrix of the 
quadratic bosonic Haniiltonian to a 48 x 48 matrix. The 
plaquette phase has been found to be stable in the do- 
main 1.64 < J < \/6 with soft modes appearing at tf = 
when J w 1.64. 



where G is a vector belonging to the reciprocal lattice of 
the lattice defined by the vectors a and b, and M{q) = 



( 



\ 



-1 








-1 



g-i9x|a| 



1 





1 





1 



. g-igx|a| 



gi9x|a| 





is the Fourier transform of the interaction matrix. The 
analysis of the second-order terms in (33) shows^'^ that 
the paramagnetic solution pi ~ becomes unstable at 

r/J = \/6 at the wavevectors qa = (eR'^)' 9b = 
(^M' 2|F|)' ^^'^ ~'1Bj indicating a transition to a 

phase of periodicity (3a — 6) x 46. 

The approach of Ref. 9 consists in keeping in the energy 
functional (33) only the critical modes with q — ±qA 
and q = ±qB whose amplitudes are described by Fourier 
coefhcients 



C. Analytical study of the critical region below Tc 



PqA 



where 



(34) 



In this subsection it will be convenient instead of 
Eq. (9) to use a different parametrization of the classical 
spins of norm S, 



s} = Sp.^ 



(31) 



In the asymptotic regime where the transverse field F 
dominates over nearest-neighbor interactions, we are in 
the polarized phase with S'f — {pi — 0). Upon decreas- 
ing the transverse field the components S'f are expected 
to deviate from zero. To sixth order in the p^'s, the clas- 
sical energy of the model is given by 



E = -jj2 m.,,p,pj-tJ2 ( 1 - t - 



16 



(32) 

Let us denote by p^^^^^ = J^qPq.ne'^'^ with n = 1, . . . , 4 
the values of pi on the four sublattices (see Fig. 6). Since 
Pj^ ^ is real, p^^ — P*-qn- energy per site £ is then 
given by 



J 



P-q,r 



M{q) 



Pq-.'. 



+ 



32 



64 



<33,94 

«, 91, 92, 93, \«=1 

94,95,96 



5 " 

9l+92+93-(-94,G 



(33) 



91 +92+93+94+95 +96, G 



UA = ( 1, e* 12 , Fe' 12 , Fe 



(35) 



are the eigenvectors of M{qA) and M{qB) associated to 
the eigenvalue \/6 and 



F 2 sin 



57r 
12 



x/2 



(36) 



In the framework of this approach, fg^\ the sum of 
the second and fourth order contributions to Eq. (33), is 
given by: 



(4) 



-\{T,-T){1 + F%\pa? + \pb\ 



-fF2 [\PA? + \PB?] 



(37) 



and depends only on \pa\'^ + \pb 

(4) 



2 7,i 



The minimum of is achieved when 



\pa\ + \pb\ - Y~ 



(38) 



from which it follows that, to leading order, 
\pa\ ~ IpbI ~ (rc - r)5 and Sl^^ ^ (Fe - Tf. However, 
condition (38) leaves both the ratio |pb|/|p^| and 
the phases ijiA and (pB completely undefined. To find 
them one has to consider also the sixth order term in 
Eq. (33),^'^ which for the critical modes reduces to 



(6) _ sr 



3r 

2 



(l-fF6) [\PA? + \PB?Y 

F3 [|p^|5|pb|cos(50a-(/)b) 
+ \pb\^\pa\ cos{5(j)B - (I)a)] 



(39) 
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The general structure of Eq. (39) has been derived in Ref. 
8 from the symmetries of the problem. 

According to the previous discussion, to leading order, 

£q^^ ~ (Fc— r)'^. For all values of the amplitudes \pa\ and 
\pb\, the expression in the right-hand side of Eq. (39) is 
minimal when both cosines are equal to —1. This selects 
the phases: 



OA = 



B — 



— V 
12^ ' 



(40) 



where p is an integer, yielding 24 independent sets 
{(j)A,4>B)- The variation of £g^^ with respect to \pa\ and 
\pb\ under the constraints (38) and (40) then selects ei- 
ther \pb\/\pa\ = F or \pb\/\pa\ - F^^. All 48 solu- 
tions thus generated correspond to the same dimer pat- 
tern (shifted or/and rotated) found in Ref. 9 and thus we 
recover the 48 fold degeneracy discussed in Ref. 8. 

The approach described above is based on the assump- 
tion that all other modes would only contribute to the 
energy expansion to higher order. We shall now show 
that, since when considering only the critical modes one 
has to push the expansion to order 6, this assumption is 
not valid because some second- and fourth-order terms 
involving noncritical modes also make contributions of 
order (Fc — F)'^ which are essential for determining 4)a 
and 4>B- 

The dominant terms coupling the critical modes with 
q = -iiqA and q = igs with extra modes are expected to 
be linear in the amplitudes of these extra modes and of 
the third order in the amplitudes of critical modes. The 
conservation of the total momentum then imposes on the 
wavevectors of these extra modes the condition: 



ruAqA + TUBqB , 



(41) 



where niA and wlb are integers and mA + "nis is odd. 
In the first Brillouin zone there are only two wavevectors 
compatible with this condition: qc = 2qA — qB and ~qc- 
Let us denote the Fourier coefficients associated to the 
modes with (f = (fc by p„ = |p„|e"^", where n = 1, . . . , 4 
refers to the number of the sublattice. The terms in the 
energy functional which are linear and harmonic in p„ 
are 



^1 



(4) 



M{qc) - -1 



-Y 

8 ^ 

n— 1 



{RnPn + C.C.) 



(42) 



Injecting Eq. (44) into Eq. (42) we obtain 



21 3 



= -Th{TlJ)[\pA\^ + \pB\^] 

-F.9 (F/J) [\pa?\pb\ cos{b(j)A - ^b) 
+ |pa||/0b|^cos(5(?!)s - 4>/ 



(45) 



where we have introduced the notation 

h{l) - 8(^7_ 3) {7(1 + F"") + 672(3F^ - 1)} , 

5(7) = ^(|^{47^^ + 3^/2(2F^-l)}. 

Eq. (44) proves that p„ scales as 

(F.~F)i, 



\PA\ 



\PB\ 



(46) 



leading to E^f' ^ (F^, — F)'^. So it is clear that this contri- 
bution cannot be neglected since it is of the same order as 
^Q^'', and that other contributions involving non-critical 
modes such as e.g. sixth order terms will be of higher 
order. This means that the phases of the critical modes 
have to be determined by minimizing the sum of £q^^ and 
S^'^K The contribution to this expression depending on 
the phases reads: 



J J 2 



[\PA\^\PB\cOs{5(j)A - (j)B) + 
+ \pa\\pb\^ C0s{5<j)B - <j)A)] ■ 



Now g{T/J) - (3/2)F3 is positive for F/J > ^3. There- 
fore, since we are interested in the domain just below 
F/J = VG, the energy is minimal when both cosines are 
equal to -|-1. This selects the phases: 



12 



12 



(47) 



with 



where p is an integer. This leads again to 24 indepen- 
dent sets {4>Aj4'b)- In addition, minimizing sjf^ + £['^^ 
with respect to the amplitudes \pa\ and \pb\ under the 
constraint (38) selects, as before, either |pB|/|/9yi| = F 
or |/9_b|/|pa| = F~^. The 48 resulting solutions corre- 
spond to the 24 equivalent dimer patterns which can be 
obtained from the one shown in Fig. 7(a). The difference 
between Eq. (40) and Eq. (47) explains the qualitative 
difference between the structures of the plaquette phase 
found in this work and the solution of Ref. 9, which does 
not disappear even when the amplitudes of the q = ±qc 
modes become negligible as compared to those of the 
critical modes. 



i?„ = pUuAt + i{p*A)Hu*A)lpB{uB)^ 

+3pAiuA)M?i^B*)l + pI{ub)1 ■ 
The variation of Eq. (42) with respect to p* gives 



F 



Pn 



2J^ 



M{qc) 



-J 



(43) 



(44) 



IV. INTERMEDIATE MIXED PHASES 

During the numerical minimization of the classical en- 
ergy for the 192-site system with periodic boundary con- 
ditions an additional intermediate phase was found to ex- 
ist between the columnar and the plaquette phases. We 
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refer to this intermediate phase as the mixed phase be- 
cause in the dimer representation the bonds with larger 
dimer densities, dij > ^, are arranged in an alternat- 
ing pattern of plaquettes and columns (see Fig. 8(a)). 
The mixed and plaquette phases have the same transla- 
tional symmetries. However, the point group symmetries 
of the gauge-invariant dimer patterns in the two phases 
are different: P31m for the plaquette phase (see Figs. 
7(a)) and Cmm for the mixed phase (see Fig. 8(a)). The 
phase transition between these two phases has to be of 
the first order, since the symmetry groups are not such 
that one is a subgroup of the other. 




(a)Mixed phase: dimer 
representation 




(b)Mixed phase: spin 
representation 

FIG. 8: (Color online) (a) The dimer density dij in the mixed 
phase at T / J = 1.72. The thickness of the bonds is propor- 
tional to dij . The dimer densities are also emphasized by the 
colors of the bonds ranging from red (the lowest densities) to 
dark blue (the highest densities, dij > |). The bonds with 
dij > I are organized in an alternating pattern of plaquettes 
and columns, (b) The spin configuration in the same state 
in the gauge that leads to a 24-site unit cell (large hexagon). 
Note the existence of 2 sites at which Sf — 0. 

As in the case of the plaquette phase, the size of a unit 
cell of the mixed phase can be reduced from 48 sites for 
the standard gauge shown in Fig. 1 to 24 sites in the 
gauge of Fig. 7(c), see Fig. 8(b). In this gauge the spin 
pattern consists of spins with the same sign of S'f hav- 
ing seven different orientations, one of which is in the 
direction of the field. In contrast to the spin pattern in 
the plaquette phase, which is centered on one of the sites 



of the honeycomb lattice (Fig. 7(c)), in the mixed phase 
this pattern is centered on one of the bonds of the lattice 
(Fig. 8(b)), which explains the difference in symmetry be- 
tween the two states. The degeneracy of the mixed phase 
is equal to 36 in terms of the dimer representation and 
to 72 in terms of the spin representation. Each of the 36 
equivalent dimer patterns corresponds to two spin config- 
urations which can be transformed one into the other by 
changing the sign of S'f for all spins. The stability of the 
mixed state with respect to small fluctuations has been 
investigated with LSWT in the gauge producing a 24-site 
unit cell, and this phase has been found to be stable in 
the range 1.394 < T/J < 1.774. 

The existence of the mixed state whose structure is 
shown in Fig. 8 suggests that there can also exist states 
in which the straight rows of plaquettes are still equidis- 
tant but separated not by single columns but by a larger 
number of columns which below is denoted by n (see 
Fig. 9). From now on we number such mixed states by 
the index n and call the simplest mixed state discussed 
above the first mixed state. 
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FIG. 9; (Color online) Dimer patterns in the mixed states 
with n < 3, where n denotes the number of columns separat- 
ing the plaquette patterns. The notation is the same as in 
Figs. 7(a) and 8(a). The bonds with dij > | are organized in 
an alternating pattern of plaquettes and columns. Thin lines 
show the boundaries between unit cells. 

It is not hard to understand that the unit cell of the 
second mixed state (in the optimal gauge in which the 
sign of is the same for all spins) has exactly the same 
symmetry as the unit cell of the first mixed state and 
can be obtained from it by adding on each side eight 
more sites. The successive repetition of this procedure 
allows one to construct the unit cell for any integer n 
and to find that it contains 8(2?! -I- 1) sites. However, 
due to the symmetry of the unit cell, the number of non- 
equivalent sites only increases by four when n increases 
by one, which leads to 4n -f 3 non-equivalent sites. 

For n < 7 we have performed a numerical minimiza- 
tion of the energy for the unit cells corresponding to such 
structures, and we have found that, upon decreasing T/J, 
the energy of the second mixed state first becomes lower 
than that of the first mixed state, after what the energy 
of the third mixed state becomes lower than that of the 
second mixed state, and so on. Tab. I summarizes the 
values of T/J at which the transition between the nth 
and (n -I- l)th mixed states takes place and reports the 
width of the region in which the nth mixed state has the 
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lowest energy. It can be seen that for n > 1 this width 
is scaled down by a factor of the order of 50 each time n 
increases by 1. This means that r"'"+^ approaches a fi- 
nite limit exponentially fast. The extrapolation shows 
that the accumulation point of r"'"+^ at n — >■ oo is 
Tf /J = 1.67612786261. Below this field columnar states 
have the lowest classical energy. 



n 




Ar"/j 


1 


1.69372479498 


4.3 X 10"'' 


2 


1.67655449242 


1.7 X 10~^ 


3 


1.67613666486 


4.2 X 10"^ 


4 


1.67612802659 


8.6 X 10"^ 


5 


1.67612786551 


1.6 X lO"'' 


6 


1.67612786267 


2.8 X lO"'' 


7 







TABLE I: Critical field r"'"+^ of the transition between the 
nth and (n+l)th mixed states. The last column shows AF" = 
pn-i,ii _ pn,7i+i^ the field range in which the nth mixed state 
has lower classical energy than the (n— l)th and the (n-|- l)th 
states (r^'^ refers to the transition between the plaquette and 
the first mixed state). 

Note that it was impossible to discover any of the mixed 
states with n > 1 during the minimization of the energy 
for the 192-site cell (with periodic boundary conditions 
and the standard gauge of Fig. 1) which was instrumental 
in discovering the n — 1 mixed state. The reason is 
very simple - the periodicity of all the states with n > 1 
is incompatible with the periodic boundary conditions 
implemented in this 192-site cell. 

The existence of such a sequence of phase transitions 
suggests that the main contribution to the energy of the 
nth mixed phase (counted off from the energy of a colum- 
nar state) is proportional to the density of linear defects 
(vertical rows of plaquettes) whose energy can be consid- 
ered as linearly dependent on F, whereas the main cor- 
rection to this energy comes from the repulsion of nearest 
defects, which decreases exponentially fast with the dis- 
tance between them. This was checked at F = F^ where 
the proper energy of a linear defect changes sign, and in- 
deed we have found that the energies of different states 
are compatible with an interaction of linear defects that 
is exponential in the distance between them. This makes 
us confident that the narrow region above T°° has to con- 
tain an infinite sequence of mixed phases with all integer 
indices n. 

It is well known that in a system consisting of a se- 
quence of linear defects there can also exist phases with 
more complex structures, in which the linear defects are 
not equidistant. In terms of our problem such phases 
would correspond to a regular alternation of, for exam- 
ple, n and n + I columns, or of n, n and n -\- 1 columns, 
etc., leading to what is known as a devil's staircase. ^'^ 



Usually such phases appear in a phase diagram if the in- 
teraction of more distant defects is also repulsive, whereas 
when the interaction between next-to-nearest defects is 
attractive, one gets a direct transition from the nth to the 
(n + l)th phase without the presence of an intermediate 
(n, n -f 1) phase. 

We have verified numerically that in our system the en- 
ergy of the (1, 2) mixed state is never lower than either 
the energy of the first state or that of the second mixed 
state, which means that it cannot be present in the phase 
diagram. Quite surprisingly, the situation with the (2, 3) 
phase is different, and in a narrow interval around F^'^ 
[from F2'3 - 1.3 x lO'^ to F^'^ + l.g x IQ-^] its energy 
is lower than the energies of the second and third mixed 
states. One can estimate that even if some other com- 
plex phases do exist, the field range where any of them 
minimizes the energy will be at least a couple of orders 
of magnitude smaller than the already extremely narrow 
interval of the existence of the (2, 3) state, so we decided 
not to pursue the investigation of this point any further 
since it cannot be of much relevance. 

A more important question is whether the plaquette 
and the first mixed states may be separated by a region 
where there appear mixed states of a different type, in 
which the density of columns is lower than in the first 
mixed state, so that the neighboring columns are sep- 
arated by domains of plaquette state. Such a scenario 
seems to us to be impossible however for the following 
reasons. 

The comparison of Fig. 7(c) with Fig. 8(b) suggests 
that the structure of the first mixed state is very close 
to what one would obtain by constructing the superpo- 
sition of two plaquette states centered on neighboring 
sites of the lattice (and letting this superposition relax) . 
Therefore one can interpret these two states as different 
manifestations of a unique state which can move around 
in a complex periodic potential with minima both at the 
positions corresponding to lattice sites and at the posi- 
tions corresponding to the middles of lattice bonds. For 
F > F^^i = 1.73690830184J, the minima located at lat- 
tice sites are the lowest, whereas for F < F°^^, the min- 
ima located at the middle of lattice bonds are the low- 
est. Exactly at F = F^^^ all these minima have equal 
depths. This picture can be confirmed by constructing a 
family of states which continuously interpolates between 
the plaquette and the first mixed state, which allows 
a numerical analysis of the effective potential discussed 
above. This analysis reveals that at F = F"'^ the barrier 
separating unequivalent (but equal) minima is very low 
(~ 1.07 X lO^^J per site). Nonetheless, any attempt to 
construct a state which somewhere looks like the plaque- 
tte state and elsewhere like the first mixed state would 
force the system to overcome this barrier in some places. 
This will increase its energy in comparison with that of 
the plaquette or of the first mixed state. 

The numerical evidence in favor of this conclusion 
comes from observing that the state which would differ 
from the first mixed state by having half its density of 
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columns has a periodicity which is compatible with the 
192-site cell used in our numerical energy minimization. 
Therefore, if at F = r°'^ the energy of this state was 
lower than that of the plaquette and of the first mixed 
states, this state would be accessible during this mini- 
mization procedure. To be on the safe side, we have also 
performed a minimization of the energy for the cell whose 
periodicity in addition to the formation of the plaquette 
and of the first mixed states allows for the appearance of 
the states which differ from the first mixed state by keep- 
ing only one column out of three (or two out of three), 
but this has not allowed us either to find any state with 
energy lower than that of the plaquette or of the first 
mixed state. This gives an additional evidence in favor 
of our conclusion that the phase transition between the 
plaquette and the first mixed states should be a direct 
one without any intermediate phases with a more com- 
plex structure. 



V. PHASE DIAGRAM 

A. Classical phase diagram 

The classical phase diagram consists of 4 regions: (i) 
The columnar phase, which is highly degenerate since 
all columnar states have the same energy. It extends up 
to r/J w 1.676; (ii) The region of mixed states with 
columnar patterns separated by straight rows of pla- 
quettes in the interval 1.676 < T/J < 1.737; (iii) The 
plaquette phase, with a 24-site unit cell, in the range 
1.737 < r/J < « 2.45; (iv) The fully polarized phase 
with all spins pointing in the direction of the field for 
r/J > \/6. The transition from the fully polarized phase 
to the plaquette phase is a second-order one, all other 
transitions being of the first order. These results are 
summarized in Fig. 10. 
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FIG. 10: (Color online) Classical phase diagram in the dimer language (above) and in the spin language (below). In the dimer 
representation the thickness of the bonds is proportional to the dimer density. Thick blue bonds correspond to the highest 
dimer density. In the spin representation the radii of the circles are proportional to Sf and arrows indicate the orientation of 
the classical spins. 



B. Quantum fluctuations 

Quantum fluctuations can a priori modify this phase 
diagram in two main ways. First of all, if the degener- 
acy of the classical ground states is accidental (that is, 
not related to symmetry) , they can select some of these 
states. This is indeed the case in the columnar phase, 
where the columnar states with domain walls of only the 
first type are selected already at the level of harmonic 
fiuctuations. 

Secondly, quantum fiuctuations can shift the phase 
boundaries. When one takes into account only the har- 
monic fluctuations, this applies only to first-order tran- 
sitions. Indeed, at a first-order transition, the classical 
energy is the same for the two competing configurations, 
but the spectra of harmonic fiuctuations are different. 



and one phase will in general be stabilized over the other 
by zero point fluctuations. A convenient way to keep 
track of the stability of the various phases with respect 
to quantum fiuctuations is to draw a phase diagram in 
the {T/J,l/S) plane (see Fig. 11) showing which phase 
has the lowest total energy. 

The resulting phase diagram can be quite involved 
when there are many phases in competition, and this 
is clearly the case here since, for 1/S — 0, there exists 
an infinite sequence of mixed phases. However, it turns 
out that for 1/S above lO^'^ only three of them survive, 
as is shown in Fig. 11. All other mixed phases exist 
only for 1/5 ^ 10~^ in a very narrow range of trans- 
verse magnetic field of width < 10~'*J. They are thus 
invisible on the scale of Fig. 11, which has been adjusted 
to properly describe the competition between the two 
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main phases (plaquette and columnar). On that scale, 
the phase diagram consists of six phases: the polarized 
phase, the plaquette phase, the first, second and fourth 
mixed states, and the columnar phase. The general trend 
is that the plaquette phase is stabilized by quantum fluc- 
tuations over the mixed phases as well as the columnar 
phase. 

Note that the transition between the plaquette and the 
columnar phases cannot be followed below V / J — 1.64 at 
this level of approximation because the plaquette phase 
is no longer locally stable with respect to harmonic fluc- 
tuations. The continuation of this boundary by a dashed 
line in Fig. 11 is just a guide to the eye. To follow this 
line further would require to go beyond the harmonic ap- 
proximation. The transition between the plaquette and 
polarized phases being of the second order, the boundary 
has to start vertically since, at the transition, both states 
have the same quantum corrections in the harmonic ap- 
proximation. This is indicated by a vertical dashed line 
in Fig. 11. To find the curvature of this line would re- 
quire to go beyond the harmonic approximation. 
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FIG. 11: Semiclassical corrections to the phase diagram, Mn 
n e {1,2,4} denote l""', 2""* and 4"" mixed states, zoommed 
for values of the field close to the accumulation point of mixed 
states. 

In view of the very strong modification of these phase 
boundaries upon decreasing S*, it is legitimate to wonder 
about the fate of the columnar and mixed phases for 
S — 1/2, for which the model can be mapped onto the 
QDM in the limit T/ J 0. The results presented above 
suggest that the mixed phases have absolutely no chance 
to extend to S = \. 

Regarding the competition between the columnar and 
the plaquette phases, we can get an estimate of the criti- 
cal value of the spin at which the boundary between them 
crosses the axis F = by looking at the linear 1/ S correc- 
tions starting from the point where the two phases have 
the same classical energy, a point that does not appear 
on the phase diagram of Fig. 11 since it lies inside the 
1*** mixed phase. This leads to the conclusion that the 
columnar phase disappears above l/S* « 0.67, i.e. below 
S K, 1.49. Note that this should probably be considered 
as a lower bound in terms of S since the boundary is 
slightly concave. So, for S = 1/2, the semiclassical cal- 



culation at the harmonic level predicts only two phases: 
a plaquette phase up to F/ J = \/&, and a polarized phase 
above. The fact that we find the point F/J = to be in 
the region of stability of the plaquette phase is in good 
agreement with the QDM, which has been found by QMC 
to be in the ^12 x ^/\2 phase at V/t = 0.^'^ 



VI. CONCLUSIONS 

In conclusion, we have investigated the classical phase 
diagram of the FFTFIM on the honeycomb lattice and 
how it is modified by the semiclassical corrections in- 
duced by harmonic fluctuations. As compared to what 
has been already known about the model, namely that 
the paramagnetic phase is unstable at F/J = ^/fj towards 
a crystalline phase with a large unit cell, the classical 
phase diagram turns out to be surprisingly rich, with 
a multitude of additional phases: a columnar phase at 
small transverse field and an infinite cascade of phases 
of mixed columnar and plaquette character. The phase 
towards which the paramagnetic phase is unstable at 
F/J = \/6 has been found to have the same symme- 
try and periodicity as the state proposed in Ref. 9, but 
a different structure. Both are characterized by a 24-site 
unit cell in the spin language, and by a 12-site cell on 
the dual lattice in the dimer language, but the state we 
have found has a plaquette structure. At the classical 
level, the columnar phase is fully degenerate, all colum- 
nar states having rigorously the same classical energy. 

Quantum fluctuations have been found to modify this 
phase diagram in two important respects. First of all, 
harmonic fluctuations have been shown to partially lift 
the degeneracy of the columnar phase in favor of the 
columnar states with only one type of domain walls. 
Since the remaining degeneracy is not related to a sym- 
metry of the model, anharmonic corrections are expected 
to lift further this degeneracy. Secondly, they modify 
strongly the phase boundaries, and for the ultra quan- 
tum limit, S — 1/2, they predict that the plaquette phase 
survives down to F — > 0. 

Going back to the original motivation of this investiga- 
tion, namely the properties of the QDM on the triangu- 
lar lattice, these results deserve a number of comments. 
First of all, our semiclassical approximation predicts that 
the phase which is the analog of the \/Vl x VT2 phase of 
the QDM has a 4-site plaquette structure. This reopens 
the issue of the nature of the y/\2 x -\/l2 phase of the 
QDM. According to the results of GFQMC simulations,^ 
possible structures are constrained by a quasi-extinction 
of the dimer density correlation function at the corner of 
the Brillouin zone. This has been shown to be consis- 
tent with a uniform distribution of dimer density inside 
the interior part of the 12-site hexagonal unit cell, a con- 
clusion somehow supported by the conclusions of Ref. 9 
regarding the nature of the phase close to the paramag- 
netic phase. Now that we know that this phase is in fact 
a plaquette phase, it would be interesting to revisit the 
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GFQMC results to see to which extent a plaquette phase 
of this type might be consistent with the quasi-extinction 
at the zone corner. 

It is also inspiring that a columnar phase appears in the 
classical solution of the FFTFIM since a similar phase is 
present in the QDM for attractive interactions between 
dimers. We did not manage to find a convincing con- 
nection between large S in the FFTFIM and negative 
V in the QDM, but since we found intermediate phases 
between the columnar phase and the plaquette phase in 
the FFTFIM, it is tempting to speculate that such phases 
may also exist in the QDM. 



states in which is minimal, that is, is equal to zero. 
In terms of dimer models such states are usually called 
staggered or nonflippable states, -"^ because they do not 
contain flippable pairs of dimers. 

Since in a staggered state all frustrated sites have iden- 
tical environments (with exactly one frustrated neighbor) 
and all non-frustrated sites also have identical environ- 
ments (with exactly two frustrated neighbors), such a 
state can be described by the same two variables 9i and 
6*2 introduced in Sec. II B for the description of a colum- 
nar state. In terms of 9i and 62 the energy of a staggered 
state can be written as 
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-E'st = — 



N 



^ (cos^ 01+4: COS 9i cos 02 - COS^ 02 
+r(sine'i +sin6l2)] . 



(Al) 

Even without minimizing Est with respect to 9i and 02 
one can notice that for any 9i and ^2 



Appendix A: Comparison of columnar and staggered 
states 

Columnar states are the states which maximize Nd, 
the number of pairs of frustrated bonds situated at the 
smallest possible distance from each other [as shown in 
Fig. 2(d)]. In this Appendix we want to compare the 
classical energy of these states with the energy of the 



Est{9u92)-E,oii9i,92) = (JiV/4)(cos^^l -COS02)' > 

(A2) 

and therefore the energy of a staggered state [the mini- 
mum of -Est (^1,^2)] has to be higher than the energy of 
a columnar state [the minimum of i?coi(^i,^2) achieved 
when cos 9i ^ cos ^2]- This proves that the maximization 
of Nd is always a better strategy than its minimization, 
even when the ratio F/J is not small. 
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